Due to the particular properties, lead monoxide (PbO) has been used as an industrial material for various technological applications and as a constituent in electronic devices, paints, gas sensors, x-ray cathodes, electrophotography, laser technology, pigments and storage batteries, radiological medical protective clothing, ointments and plasters, etc.

Two existing polymorphs of PbO with respect to temperature and pressure are α-PbO (litharge) and β-PbO (massicot).^[1–5] The red α-PbO with a tetragonal two-dimensional layer structure is stable at ambient pressures and lower temperatures while the yellow β-PbO has an orthorhombic chain structure and is stable at higher temperatures. The temperature around 860 K is the transition temperature between two phases of PbO. This transition temperature has been studied by many researchers and there are some discrepancies about this value which originates from slow kinetics of the reaction. The stability of these phases is susceptible to impurities and temperature.

According to other studies,^[6] the simulation and modeling methods are considered as supplements to experimental studies, and can be used instead of traditional time consuming and trial and error methods. These methods can help our understanding of various properties of materials under specific conditions particularly high temperatures and pressures.

The evolutions of thermodynamic and structural properties at high pressures and high temperatures are of significant importance in modern technologies which experimentally encounter with difficulties, since it is difficult to measure the mentioned properties under these conditions. Furthermore, the change in volume due to changes in temperature and pressure, eventuates in phonon property variations which relate to elastic constants and some properties of crystals such as entropy, specific heat, etc. High pressure results in the phase transition and the variation of mechanical and thermodynamic properties of solids. Further understanding of bulk modulus, Poisson’s ratio, Young’s modulus, and elastic constants is important for finding mechanical properties. In order to shed more light on experimental work, many theoretical studies by various research groups have been accomplished. Most of them used ab initio studies based on the traditional density functional theory (DFT) for investigating structural and thermodynamic properties of materials in the wide range of pressures and temperatures.^[7–14]

In order to investigate the effects of pressure and temperature on both phases of PbO and the phase transition between them, some theoretical and experimental studies for lead monoxide have been performed at various temperatures and pressures. Few theoretical ones are performed on the properties at high temperature and pressure using first principles phonon calculations. The methods to investigate these two systems are still a topic of debate. Boher et al.^[15] reported a transition from alpha phase to an intermediate phase with an orthorhombic Cmma structure named γ-phase which is stable below 200 K. In this phase the structure is a little distorted. Whereas the difference in the atomic site coordinates between these two structures is almost 0.002 nm, so in most studies this phase is ignored, but according to Adams et al.,^[16] α-PbO phase transition to the mentioned distorted phase occurs under pressures 0.7 GPa and 2.5 GPa. This structure is an intermediate phase between alpha and beta phases of PbO. The alpha–beta phase transformation is slow and pressure sensitive. According to White et al.,^[17,18] these two phases are in equilibrium at 1 bar (1 bar = 1.0 × 10⁵ Pa) and 813.15 K. In both α-PbO and β-PbO phases, the Pb²⁺ cations are close to one another provided that anions are placed at various locations. Studies proved that these two phases possess a layer structure and some specifications of molecular crystals due to inert pairs of cations. Canepa et al.^[19] have investigated the electronic, structural and dynamic properties of these two phases using LCAO approach within density functional theory framework with various GGA and hybrid functional approximations. They introduced dispersive interactions for calculating cell parameters as a major component of the interaction between layers and demonstrated that the anisotropy of Pb lone-pair orientation within both phases, affects α-PbO more than β-PbO.

In this study, we focus on the structural and thermal properties of β-PbO under high pressure and high temperature employing density functional theory and quasi harmonic approximation. All the calculations have been performed using Debye and the Debye–Grüneisen model.^[20] Bulk modulus, specific heat, thermal Debye temperature, thermal expansion coefficient, Grüneisen parameter of beta phase of PbO have been compared.

We have employed density functional theory (DFT) as implemented in the Quantum Espresso package^[21] incorporated with plane-wave basis set and ultrasoft pseudopotentials^[22] to study the temperature and pressure dependence of thermodynamic properties related to β phase of lead oxide. Exchange and correlation effects are treated by the GGA approach with Perdew–Burke–Eruzerhof (PBE) functional^[23] and core electrons are represented explicitly by using Vanderbilt ultrasoft pseudopotentials with a plane-wave function cutoff and energy cutoff of 75 Ry (1 Ry = 13.6056923(12) eV) and 600 Ry, respectively. In order to confirm the convergence of calculations, we attentively verify k-points set, so a grid of 12 × 12 × 12 Monkhorst–Pack^[24] of k-points was considered. Besides, the ionic positions and cell shapes of β-PbO are relaxed by employing the Gaussian smearing factor of 0.05 eV. For phonon calculations 12 × 12 × 12 k mesh in the BZ and 80 Ry energy cutoff have been used.

To investigate the thermodynamic properties and phase stability of solids at high temperature and high pressure, we applied a very quick and simple model, the quasi-harmonic Debye model in which the vibrations in a crystal behaved as a gas of non-interacting phonons. For a solid which is described by the relationship of energy and volume, the Gibbs program can estimate the Debye temperature Θ(V) in order to calculate Gibbs free energy and equation of state (EOS) and other thermodynamics relations. The Birch–Murnaghan family of equation of state (EOS) which comes from assuming a polynomial form for the energy, is used for illustration of pressure–volume curves by a third order,^[25–27] defining, B = −V(dp/dV); B′ = (dB/dp) and x = (V/V), the third order Birch–Murnaghan takes the form 1 2 where B₀ is bulk modulus at zero pressure.

According to standard thermodynamics the non-equilibrium Gibbs free energy G*(V; P, T) of the crystal phase can be written as follows: 3 where E(V) is the static energy and F_Vib is the non-equilibrium vibrational Helmholtz free energy.^[28–33] Minimizing Eq. (3) with respect to volume gives the following mechanical equilibrium condition: 4 where p_sta is the static pressure, p is the applied external pressure and p_th is thermal pressure as follows: 5 6 Using Debye model of phonon density of states, the vibrational term F_Vib which includes the vibrational contribution to the internal energy and the −TS constant temperature condition term is given as 7 8 where the vibrational frequencies of phonon modes ω_j depend on the crystal geometry. The Helmholtz free energy versus volume has been fitted to isothermal third-order finite strain equations of state.

In quasi-harmonic Debye model, the adiabatic bulk modulus is equal to the isothermal bulk modulus B_T. Thus it can written in the following equation: 9 where E is the total energy of crystal at 0 K and V_opt(p,T) are the volumes which minimize the non-equilibrium Gibbs function G*(V; P,T) that gives the thermal equation of state (EOS) and the chemical potential μ (p,T).

When the volumes at equilibrium are known, some thermodynamic properties can be derived from Eq. (8). Considering this and the isothermal bulk modulus B_T, Helmholtz free energy (F), Gibbs free energy (G), internal energy (U), equilibrium entropy (S), heat capacity at constant volume (C_v) can be given by^[34] 10 11 12 13 14 The thermodynamic Grüneisen parameter (γ_th) is known as 15 Then the thermal expansion coefficient (α), constant-heat capacity (C_p) and adiabatic bulk modulus (B_S) are as follows: 16 17 18 With the purpose of calculating thermal properties of lead oxide, we employed Debye and Debye–Grüneisen approaches implemented in the GIBBS code.

The Debye–Grüneisen model is applicable in various cases in many senses. In fact, this model is considered as a modified Debye model, which supposes that Debye temperature is a function of density, θ_D = θ_D(ρ), and the Grüneisen parameter can be defined as 19 In other words, Debye temperature is calculated at the static equilibrium and the quasi-harmonic approximation is applied to calculate the volume dependence of the Debye frequency. Hence, the Grüneisen ratio are also chosen to be in the following form: 20 when a = −1/2, b = 1/2, one obtains the Dugdale–McDonald (DM) approximation. For thermal property investigations the mentioned models have been used with DM approaches.

The quality of energy and vibrational data used for finding thermodynamics properties ascertain the accuracy of calculated results. Since DFT method is usual for these kinds of calculations, so the source of uncertainty relates to exchange–correlation functional parts which are systematic. One of the suggestions for improving the calculations is applying corrections to the static E(V) curve and shifting it. These corrections are called empirical energy corrections (EEC) and are divided in three different corrections: the PSHIFT EEC, the APBAF and the BPSCAL EEC.^[34] In all three EECs, the values of free parameters α, Δ p, B_exp, and V_exp are chosen so that the experimental equilibrium volume at room temperature is reproduced. In this work, the PSHIFT EEC is applied to the static energy implemented in all models which is defined as 21

Using the first principles method, we made the total energy electronic structure calculations to obtain lattice parameters at equilibrium for the lead (II) oxide. Lead monoxide (PbO) exists in one of two lattice formats, tetragonal and orthorhombic. The tetragonal crystal structure is called litharge or alpha lead monoxide. Massicot or beta lead oxide is the orthorhombic variety of PbO. These two crystal structures can be changed to the other form by controlled heating and cooling. For β-phase of lead oxide, the unit cell of the corresponding crystal was considered. The β-phase lead oxide has Pbcm space group with a = 5.49 Å, b = 5.89 Å, and c = 4.75 Å.^[24]

All calculations have been shown for both Debye and Debye–Grüneisen models for comparison. In all figures “debye” and “dg” indicate calculations related to “Debye” and “Debye–Grüneisen” models, respectively. In order to determine the equilibrium structure parameters, several volumes have been varied and a set of volumes versus energies were obtained. The calculated (V, E) data were fitted to the Birch-Murnaghan thrid-order EOS for both models and for beta phase. The volume variation versus pressure at different temperatures for β-PbO is presented in Fig. 1.

Fig. 1.

	Figure Option ViewDownloadNew Window
	Fig. 1. (color online) Variation of volume V (1 Bohr = 5.29 × 10−11m) with pressure (GPa) for β-PbO.

Through the quasi-harmonic Debye model, the thermodynamic values of various temperatures and pressures are calculated using E–V data at zero pressure and temperature. The thermal properties are calculated in the temperature and pressure ranges 0–1000 K and 0–100 GPa, respectively. Calculation results are presented in the following. According to Fig. 1 and Table 1, it is concluded that there is good agreement between the calculated results and other experimental and theoretical studies which indicates that this method is suitable for studying high-temperature and high-pressure behavior of beta phase of PbO.^[30] This offers the accuracy and reliability to our deeper investigation.

Table 1.

Table 1.

Table 1. Lattice constants, volume (Bohr3), bulk modulus (B), thermal expansion (α), and entropy for β-PbO at ambient temperatures. . This study Other studies a/Å 5.49 5.489a,b b/Å 5.89 5.895a, 5.891b c/Å 4.75 4.756a, 4.775b V 1037.82 1038.3a, 1041.81b B/GPa 33.1 22.72 ± 6.2d α/K−1 4.08 × 10−5 6.86 × 10−5b S/J⋅mol−1 · K−1 78 68.7f,g a Ref. [43], b Ref. [39], c Ref. [45], d Ref. [16], e Ref. [46], [47], f Ref. [37], g Ref. [38].

Table 1. Lattice constants, volume (Bohr3), bulk modulus (B), thermal expansion (α), and entropy for β-PbO at ambient temperatures. .

The bulk modulus B relates to interatomic potentials and can be found by the second derivative of the internal energy.^[30] The pressure dependences of the isothermal bulk modulus for β-PbO for both studied models are shown in Fig. 2. According to this figure, it is concluded that increasing pressure results in augmentation of bulk modulus, at a defined temperature. For both Debye and Debye–Grüneisen models, this trend is the same. It shows that the effect of pressure on bulk modulus is very important. According to this figure, at low temperatures (0 and 250 K), bulk modulus variation with temperature is more significant in comparison with high temperatures.

Fig. 2.

	Figure Option ViewDownloadNew Window
	Fig. 2. (color online) Pressure dependence of the isothermal bulk modulus of β-PbO at various temperatures.

The temperature dependence of bulk modulus reflects the anharmonic interactions partially, provided that for a purely harmonic crystal, the bulk modulus would be independent of temperature. The calculated bulk moduli at ambient temperature are listed in Table 1, where comparison is done with other studies.

Heat capacity C_v and C_p of a substance is another thermodynamic parameter which is important and provides necessary comprehension of vibrational properties. Pressure dependence of heat capacity at constant volume C_v and at constant pressure C_p at different temperatures (T = 0, 250, 510, 760, and 1000 K) are plotted in Figs. 3(a) and 3(b), respectively, for both models.

Fig. 3.

	Figure Option ViewDownloadNew Window
	Fig. 3. (color online) (a) Heat capacity at constant volume (Cv) and (b) heat capacity at constant pressure (Cp) versus pressure for different temperatures for β-PbO.

By considering Fig. 3(a), it is deduced that C_v values decrease smoothly as pressure increases. The important issue is that for low temperatures (0 K and 250 K), pressure effect is more significant than high temperatures. This means that at 0 and 250 K, as pressure increases, C_v variation range is around 15 J/mol K, while for high temperatures (510 K, 760 K, and 1000 K), the pressure increasing does not make any significant changes in C_v values and for all these temperatures C_v is identical. As it is obvious in this figure, for pressures ranging from to 100 GPa, at ambient temperature, C_v varies from 48 J⋅mol⁻¹·K⁻¹ at 0 GPa to the value of 33 J⋅mol⁻¹·K⁻¹ at 100 GPa. In other words, it can be concluded that for temperatures higher than ambient temperature, pressure enhancement will not affect heat capacity at constant volume (C_v). For this structure C_v values for higher temperatures are around 50 J⋅mol⁻¹·K⁻¹. This behavior is the same for both studied models Debye and Debye–Grüneisen. The heat capacity at constant pressure, C_p, is presented in Fig. 3(b). It is clear that, for all temperatures, generally, C_p decreases with pressure increasing. The decreasing trend is fast till 10 GPa. For pressures higher than 10 GPa, the inclination is smoother. Considering these calculations, C_p value at ambient conditions is almost 48 J⋅mol⁻¹·K⁻¹. The theoretical value for this structure is about 45.8 J⋅mol⁻¹·K⁻¹. As a result, the discrepancy is only about 4%. However, C_p exhibits apparently different features at high temperature.

It is clear that the effect of pressure increment at higher temperatures is more considerable for C_p values in comparison with C_v values. This behavior is the same for both studied models. Furthermore, for both Debye and Debye–Grüneisen model this slope is sharp at high temperatures during pressure increasing. To the best of our knowledge, no theoretical or experimental studies are available for C_p values at high pressures for β-PbO. However, upon pressure increasing, in agreement with other experimental studies,^[35] the specific heat falls off.

Figure 4 exhibits the dependence of C_p and C_v on temperature for β-PbO. In agreement with other studies^[36] studying the behavior of heat capacity C_v and C_p for low temperatures demonstrates that, as temperature increases, C_v and C_p also increase. Provided that for the above ambient temperature, the increasing speed of C_v becomes slow and approximates to a constant value while increasing speed of C_p is more but finally also converges to a constant value. It is obvious that at high temperatures, T > 400 K, C_v does not possess much dependency on temperature and tends to the Dulong–Petit limit^[30] and is equal to 49.7 J⋅mol⁻¹·K⁻¹ in agreement with other studies.^[36,37] For intermediate temperatures, atom vibrations control C_v value. For T < 400 K, C_v depends on temperature and grows as temperature increases. By comparison of C_v values as function of temperature, similar behaviors are seen for both Debye–Grüneisen and Debye model. However, for C_p values at high temperatures, these models are not alike. In high temperatures, anharmonic behaviors cannot be ignored for these models. Previously, Spencer and Spicer^[33] reported that the C_p (cal/deg mol) of β-PbO for temperatures 298–923 K can be described by 22 According to the mentioned relation and other experimental studies^[37,38] for high temperatures, at T = 1000 K, C_p is approximately around 56.4 J⋅mol⁻¹·K⁻¹ and 57.76 J⋅mol⁻¹·K⁻¹, respectively. As it is obvious in Fig. 4, the calculated value for C_p in this study, using Debye–Grüneisen model is 55.5 J⋅mol⁻¹·K⁻¹ which is consistent with the mentioned values, whereas, Debye model result for C_p is a little more than the reported values in reference tables.^[38]

Fig. 4.

	Figure Option ViewDownloadNew Window
	Fig. 4. (color online) Heat capacity at constant volume (Cv) and heat capacity at constant pressure (Cp) versus temperature for β-PbO.

The entropy (S) variation with respect to temperature for β-PbO is displayed in Fig. 5. It is seen that with the temperature enhancement, the entropy increases which is a general trait for every material. The reported entropy for this structure at temperature 298.15 K and 1000 K are 73.37 J⋅mol⁻¹·K⁻¹ and 152.03 J⋅mol⁻¹·K⁻¹, respectively.^[38]

Fig. 5.

	Figure Option ViewDownloadNew Window
	Fig. 5. (color online) The variation of the entropy (S) with temperature of β-PbO.

According to our calculations and Table 1, at the mentioned temperatures, the entropy for β-PbO are 78 J⋅mol⁻¹·K⁻¹ and 143 J⋅mol⁻¹·K⁻¹, the discrepancies between our calculations and the experiments are 6.3% and 5.9%, respectively. It is observed that there is good agreement between these calculations and values estimated from other theoretical and experimental studies.^[38,39] For both models the entropy variation feature is the same.

Figure 6 displays the thermal expansion coefficient versus temperature at zero pressure for β-PbO using both models. We can perceive that as temperature increases, thermal expansion coefficient increases exponentially for T < 400 K, it gradually tends to linear increase for T > 400 K and the trend of increasing becomes temperate at high temperatures. At standard condition, the value of the thermal expansion coefficient for the β-PbO is about 4.08 × 10⁻⁵ K⁻¹ in agreement with these of experimental studies.^[40] This figure indicates that Debye–Grüneisen model is more reliable to study the expansion coefficient of this material in comparison with Debye model. It is obvious that with the temperature increasing, no linear behavior of the thermal coefficient is anticipated due to anharmonicity.

Fig. 6.

	Figure Option ViewDownloadNew Window
	Fig. 6. (color online) The thermal expansion versus temperature for β-PbO.

The vibration alternations of a crystal lattice can be explained by the Grüneisen parameter (γ) which is based on the change in volume, size and dynamics of lattice as a result of the temperature change. The Grüneisen parameter for β-PbO as a function of temperature at zero pressure and versus pressure for various temperatures is displayed in Fig. 7. As it is clear in Fig. 7(a), increasing temperature eventuates in Grüneisen parameter promoting uniformly. The temperature dependency of this parameter at low pressures relates to quasiharmonic approximation. However, Grüneisen parameter decreases during the pressure increasing as it is shown at various temperatures in Fig. 7(b). Note that at higher pressures, this parameter is approximately independent of temperature. To our knowledge, there is no experimental or theoretical study about the thermodynamic properties at high pressure and high temperature for β-PbO. This may be related to the lone-pair bonding situation in PbO.

Fig. 7.

	Figure Option ViewDownloadNew Window
	Fig. 7. (color online) Trends in Grüneisen parameter as a function of (a) temperature and (b) pressure for β-PbO.

The phase transition from alpha to beta phase which occurs at about 763–837 K^[41–44] is critical for inspecting the applicability of the method used in this study in order to understand the behavior of α-PbO. All the same calculations can also be done for α-PbO in the next studies using other models. The consistency between our calculations and other studies gives a good idea about the accuracy and reliability of our present calculations for more thermodynamic and elastic properties of phases of PbO.

First principles calculations within quasi-harmonic approximation have been performed on structural and thermodynamics properties of β-PbO structure over a wide range of pressures and temperatures. Two different models, Debye and Debye–Grüneisen were employed for studying thermodynamic properties. Calculated results for volume, bulk modulus, heat capacity at constant volume and constant pressure and entropy are excellently consistent with available experimental results. The results of Debye model for heat capacity at constant pressure for β-PbO indicate that this method is more reliable for calculating heat capacity.

The discrepancy between some calculated results like heat capacity and entropy with experimental values might rise from the approximation which is made at high temperature in quasi-harmonic Debye model. According to the high percentage agreement of our calculations and other studies, this method can be used for finding Gr¨uneisen parameter and other thermodynamic properties at high pressure and temperature of Massicot with reasonable accuracy. Since most of the studies for thermodynamics properties of PbO are old and few, further investigations are required for comprehending these properties under various conditions. Furthermore, these results motivate more experimental and theoretical studies on both phases of PbO, α-PbO and β-PbO. In addition, more endeavors are necessary for making the simulation results of the high pressure and temperature better.